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Quantum tunneling, a phenomenon in which a quantum state traverses energy barriers above 
the energy of the state itself, has been hypothesized as an advantageous physical resource for op¬ 
timization. Here we show that multiqubit tunneling plays a computational role in a currently 
available, albeit noisy, programmable quantum annealer. We develop a non-perturbative theory 
of open quantum dynamics under realistic noise characteristics predicting the rate of many-body 
dissipative quantum tunneling. We devise a computational primitive with 16 qubits where quan¬ 
tum evolutions enable tunneling to the global minimum while the corresponding classical paths are 
trapped in a false minimum. Furthermore, we experimentally demonstrate that quantum tunneling 
can outperform thermal hopping along classical paths for problems with up to 200 qubits containing 
the computational primitive. Our results indicate that many-body quantum phenomena could be 
used for finding better solutions to hard optimization problems. 


Quantum tunneling was discovered in the late 1920s 
to explain radioactive decay and field electron emis¬ 
sion in vacuum tubes. Tunneling plays a key role in 
Josephson junctions, molecular nanomagnets, and in 
charge and energy transport in biological and chem¬ 
ical processes [T]. It is at the core of many essential 
technological innovations such as flash memories and 
the scanning tunneling microscope. Quantum tun¬ 
neling has also been hypothesized as an advantageous 
mechanism for quantum optimization, in particular in 
systems with thin but high energy barriers [Ml- In 
classical cooling optimization algorithms such as sim¬ 
ulated annealing the initial temperature must be high 
in order to overcome tall energy barriers. As the al¬ 
gorithm progresses the temperature is gradually low¬ 
ered to distinguish between local minima with small 
energy differences. This causes the stochastic process 
to freeze once the thermal energy is lower than the 
height of the barriers surrounding the state. In con¬ 
trast, quantum tunneling transitions are still present 
even at zero temperature. Therefore, for some energy 
landscapes, one might expect that quantum dynam¬ 
ical evolutions can converge to the global minimum 
faster than the corresponding classical cooling process. 

Quantum annealing is a technique inspired by 
classical simulated annealing which aims to take ad¬ 
vantage of quantum tunneling. The performance of 
chips designed to implement quantum annealing us¬ 
ing superconducting electronics has been studied in a 
number of recent works [EH23] • Here we consider chips 
manufactured by D-Wave Systems, described in detail 
in [24] . The qubits are subject to complex interactions 
with the environment. We show that even under such 
conditions the device performance benefits from mul¬ 
tiqubit tunneling. We consider a computational prim¬ 
itive, the simplest non-convex optimization problem 
consisting of just one global minimum and one false 


(local) minimum. Quantum evolutions enable tun¬ 
neling to the global minimum while the correspond¬ 
ing classical paths are trapped in the false minimum. 
A detailed multiqubit master equation accurately de¬ 
scribes the experimental data from the D-Wave Two 
processor at NASA Ames. We study the temperature 
dependence of the probability of success for our com¬ 
putational problem. Consistent with our quantum 
models, we experimentally determine that the tem¬ 
perature dependence of the success probability of the 
D-Wave chip is opposite to the temperature depen¬ 
dence predicted by models based on classical paths 
with thermal hopping. 

The goal of quantum annealing is to find low energy 
states of a “problem Hamiltonian” 

= a) 

/i lIV 

where the Pauli matrices crjj correspond to spin vari¬ 
ables with values {±1}- The local fields { h M } and 
couplings {J^v} define the problem instance. Quan¬ 
tum annealing is charaterized by evolution under the 
Hamiltonian 

H 0 (s) = A(s)H D + B(s)H P , (2) 

where Hp = — ^ < 7 *. The annealing parameter s 
slowly increases from 0 to 1 throughout the annealing 
time t qa . Initially A(0) B( 0). With increasing s, 

A(s) monotonically decreases to 0 for s = 1, whereas 
B(s) increases. 

The problem Hamiltonian encoding the computa¬ 
tional primitive with one global minimum and one 
false minimum is depicted in Fig. [J It consists of two 
qubit cells, left and right, each with n = 8 qubits. The 
local fields 0 < /il < 0.5 and /ir, = — 1 are equal for all 
the spins within each cell, and all the couplings J = 1 
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J = 1 



FIG. 1. Graph of the problem Hamiltonian with 16 qubits 
coupled ferromagnetically with J=1 (lines). The applied 
fields are 0 < /il < J/2 (/ir = —1) for the left (right) qubit 
cell. The symmetry and strong intra-cell ferromagnetic 
coupling makes each 8-qubit cluster evolve together. 


are ferromagnetic. The spins within each cell tend to 
move together as clusters due to symmetry and the 
strong intra-cell ferromagnetic coupling energy. We 
choose /iR > |/il|, so that in the low energy states 
of Hp the right cluster is pointing along its own lo¬ 
cal field as seen in Fig. [J The difference in energy of 
the states with opposite polarization in the left clus¬ 
ter is n(J — 2/il). Choosing h l < J/2 = 0.5, the 
global minimum corresponds to both clusters having 
the same orientation, while in the false minimum they 
have opposite orientations. 

Our aim is to distinguish quantum tunneling from 
thermal activation along classical paths of product 
states (which preclude multiqubit tunneling). We 
now explain why a classical path continuously con¬ 
nects the initial global minimum to the final false 
minimum. At the beginning of the annealing process 
B(s)/A(s) < 1, and we have (a*) ~ hk B(s)/A(s) 
(because the coupling terms are quadratic in the z- 
polarizations (erf)). As h l and /ir have opposite signs, 
so will the z-projections of spins in the two clusters 
early in the evolution. To escape this path classically 
all spins in the left cluster must flip sign, which re¬ 
quires traversing an energy barrier. The barrier peak 
corresponds to zero total z-polarization of the left clus¬ 
ter. Therefore, the barrier grows with the ferromag¬ 
netic energy of the cluster (ro/2) 2 J. The barrier height 
is much greater than the residual energy which grows 
with n(J — 2/ir)- 

In order to give a more precise description of the 
classical paths of product states, let each qubit be 
represented by a spin vector in the xz- plane. De¬ 
note by 9p the angle of the spin vector for qubit p 
with the x quantization axis. We can gain an intu¬ 
itive understanding of the effective energy landscape 
if we assume that all the qubits in the left (right) clus¬ 
ter have the same angle $l (0r). This assumption is 
based on symmetry and the strong intra-cluster ferro¬ 
magnetic energy. The resulting energy potential can 
be derived using more formal methods, like the Villain 
representation [25]. Figure [2] plots the effective energy 
potential for the left cluster as a function of 0 r with 
/il = 0.44. The classical path (red line) which follows 
the local minimum of this effective energy potential 
gets trapped in a false minimum and fails to solve the 
corresponding optimization problem, as explained in 
the previous paragraph. 

In the absence of quantum tunneling, the global 
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FIG. 2. Energy potential using /il = 0.44, plotted ver¬ 
sus annealing s and tilt angle #l of each spin vector in 
the left cluster. The red line corresponds to a path that 
starts in the initial global minimum and follows the instan¬ 
taneous local energy minimum. A second local minimum 
(dashed blue line) forms at the bifurcation point s = 0.18. 
The global minimum is found in this second path after 
s = 0.24 (dashed to continuous blue line). To reach this 
global minimum the system state has to traverse the en¬ 
ergy barrier between them (dashed green line), either by 
thermal activation or by quantum tunneling. 


minimum could be reached through thermal excita¬ 
tions along classical paths for over-the-barrier escape 
from the false minimum. This thermal activation re¬ 
sults in an increasing probability of success with ris¬ 
ing temperature. This intuition is supported by spin 
vector Monte Carlo (SVMC), a numerical algorithm 
consisting in thermal Metropolis updates of the spin 
vectors [2D]. Figure 5a confirms the thermal activa¬ 
tion in SVMC. This is opposite to both open quantum 
system theory and experiments with the D-Wave chip, 
which show a reduction of the probability of success 
with rising temperature, as explained later. Further¬ 
more, Figure [5b| shows that the probability of success 
for SVMC is lower than the probability of success for 
D-Wave and open system quantum models. 


Quantum mechanically, the system evolution goes 
through an “avoided-crossing” where the two low¬ 
est eigenstates Ei(s) and E 0 (s) approach closely to, 
and then repel from, each other (see inset in Fig. |3j) . 
Higher energy states remain well separated during the 
evolution. This level repulsion occurs due to the col¬ 
lective tunneling of qubits in the left cluster between 
the opposite ^-polarizations. At the point where the 
gap M2io(s)=-Ei(s)-.Eo(s) reaches its minimum the 
corresponding adiabatic eigenstates are formed by the 
symmetric and anti-symmetric superpositions of the 
cluster orientations. The size of the minimum gap is 
varies with h \ J , as seen in Fig. [3] 

Under realistic conditions, a quantum annealer can 
be strongly influenced by coupling to the environ¬ 
ment, for which we introduce a detailed phenomeno¬ 
logical open quantum system model. We shall assume 
that each flux qubit is coupled to its own environ¬ 
ment with an independent noise source; this is con- 
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FIG. 3. Inset shows the quantum energy gap = 

Ei(s) — E 0 {s) versus s, using /il = 0.44. The dashed 
line is the gap in the diabatic (pointer) basis. In the main 
plot, the minimum gap decreases with h l. The horizontal 
boundary of the red-filled area (324 MHz) corresponds to 
15.5 mK, the lowest temperature in our experiments. The 
lower inset shows the spin configurations of the two lowest 
eigenstates at the end of the annealing. 


sistent with experimental data |10j . The coupling of 
the environment to each flux qubit is through flux 
fluctuations, and is proportional to a a z qubit oper¬ 
ator. The properties of the noise are determined by 
the noise spectral density S(ui), which is character¬ 
ized by single-qubit macroscopic resonant tunneling 
(MRT) experiments in a broad range of biases (0.4 
MHz — 4 GHz) and temperatures (21 mK — 38 mK) 
for tunneling amplitudes of a single flux qubit below 
1 MHz . The MRT data collected is surprisingly well- 
described El EZ] by a phenomenological “hybrid” 
thermal noise model S(u>) = Rif(w) S 0 h (<*>)■ Here 
5 0 h(w) = h 2 'quie~ u,Ta /{l- e ~ hui / kBT ) denotes the high- 
frequency part, and has Ohmic form with dimension¬ 
less coupling r] and cutoff frequency l/r c (assumed 
to be very large). The low-frequency part Rif is of 
the 1 // type ESI and in current D-Wave chips this 
noise is coupled to the flux qubit relatively strongly. 
Its effect can be described with only two parame¬ 
ters: the width W and the Stokes shift e p of the 
MRT line El • The experimental shift value is related 
to the width by the fluctuation-dissipation theorem 
(e p = HW 2 /2ksT ) and represents the reorganization 
energy of the environment. The values of the noise pa¬ 
rameters measured at the end of the annealing (s = 1) 
for the D-Wave Two chip are W/(2it) = 0.40(1) GHz 
and r) = 0.24(3). 

In the analysis of the transitions between the states 
we start from the initial (gapped) stage when the 
instantaneous energy gap M2io(s) between the two 
lowest eigenstates |^>o(s))> |^i( s )) is sufficiently large 
compared to the linewidth HW. Then the coupling to 
the environment can be treated as a perturbation and 


the transition rate between these states is then given 
by Fermi’s golden rule Ti_ > o(s) ~ ®(s)S(Qio(s))/h 2 . 
Here 

2 n 

°( s ) = \(Ms)W P \Ms))\ 2 (3) 

M=1 


is a sum of (squared) transition matrix elements be¬ 
tween the two eigenstates. 
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FIG. 4. Solid lines correspond to the modeled popula¬ 
tion of the lowest energy eigenstate along the quantum 
annealing process using h,L = 0.44 at 15.5 mK (top line) 
and 35 mK (bottom line). Dashed lines correspond to 
the thermal equilibrium population. In the thermalization 
phase (red) the transition rate is fast and the population 
remains close to thermal equilibrium. As the multiqubit 
energy barrier increases, the transition rates are exponen¬ 
tially reduced with s, as shown in the inset. We define the 
slowdown regime (blue) as f 9a ri_>o < 10, and the frozen 
regime (gray) as t^Fi-m < 0.1. Comparing the data at 
15.5 and 35 mK, we see a small change in the transition 
rate relative to the larger change in the thermal equilib¬ 
rium ground state population. Therefore, the probability 
of success is lower at higher temperature. 

In the minimum gap region the (squared) matrix 
element <d(s) for the transition rate is large, and the 
system is thermalized (see Fig. |4j) . More precisely, we 
have Tf-m 1 /t qa , where the inverse of the annealing 
time 1 /t qa is an approximation for the annealing rate. 
The ground state population is given by the Boltz¬ 
mann distribution at the experimental temperature. 

After the avoided-crossing region we observe a steep 
exponential fall-off of the matrix element o(s) with 
s, eventually causing multiqubit freezing (see Fig. |4|. 
Multiqubit freezing is quite distinct from single qubit 
freezing. Single qubit tunneling m decays slowly as 
the magnitude of the transverse field H(s) decreases. 
The multiqubit transition rate, however, decays expo¬ 
nentially fast (see inset of Fig. [4]). This is due to the 
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increasing effective barrier width (see Fig. |2j) , which 
results in an exponential decrease of quantum tunnel¬ 
ing and in a slowdown of the transition rate Id^o- 
Formally, the barrier width corresponds to the Ham¬ 
ming distance 


2 n 

K s ) = H KV’oKIV’o) - (^il^l^i)l 2 / 4 ( 4 ) 

M= 1 

between the opposite ^-orientations of the left cluster 
in the two lowest energy eigenstates. The exponen¬ 
tial sensitivity of multiqubit tunneling to the width 
or Hamming distance lh(s) is the cause of the expo¬ 
nential decay of the matrix element o(s), and of the 
multiqubit freezing. 

We distinguish a slowdown phase (roughly 0.1 < 
tgoTi-^o < 10) and a frozen phase (t^aTi-m < 0.1). 
In the frozen phase, there are no dynamics. Part of 


the system population remains trapped in the excited 
state |'i/’i(s)) corresponding to the false minimum of 
the effective potential until the end of the quantum 
annealing process (see Fig. [4]). 

The success probability of quantum annealing is 
(roughly) determined by the thermal equilibrium 
ground state population during the slowdown phase. 
When the temperature grows, the ground state popu¬ 
lation decreases appreciably, while the transition rate 
changes little (see Fig. [4]). This results in the observed 
thermal reduction (see Fig. 5aI. 

When the energy gap is similar to (or smaller than) 
the noise linewidth W the environment cannot be 
treated as a perturbation. We develop a multiqubit 
non-perturbative analysis in the spirit of the Non¬ 
interacting Blip Approximation (NIBA) [3D] that cov¬ 
ers all QA stages. In the slowdown phase, when 
the Hamming distance approaches its maximum value 
Ih ~ n, the instantaneous decay rate of the first excited 
state takes the form 


ri_M) — 


dr e 


iQ. w r—h(ie p T+(Wr) 2 /2) 


7TT C (r - ir c ) 

~W csch ^J^~ 


D(t), 


(5) 


where l/r c is the Ohmic noise cutoff frequency, /3 is 
h/ksT, and the factor D(t) is provided in [24]. The 
dependence on the annealing parameter s is implicit. 
The factor D{r) is related to the tunneling perme¬ 
ability of the potential barrier in Fig. [2] (similar to the 
coefficient a). The above expression describes collec¬ 
tive tunneling of the left qubit cluster assisted by the 
environment. The crucial difference from the single 
qubit MRT theory [27, 29 is that the parameters of 
the environment in the transition rate are rescaled by 
the barrier width or Hamming distance lh(s). The ef¬ 
fective low-frequency noise linewidth is lh 1/,2 (s)W(s), 
the reconfiguration energy is h(s)e p (s) and the Ohmic 
coefficient is th(s)? 7 (s). This is important at the late 
stages of quantum annealing when Ih ~ n 1. 

We observe a very close correspondence between the 
results of the analysis with the NIBA Quantum Mas¬ 
ter Equation for the dressed cluster states and the D- 
Wave Two data displayed in Fig. [5b] We emphasize 
that for NIBA (and the standard Redfield equation 
with S’oh (a;)) we do not have any parameter fitting: 
the parameters are obtained from MRT experiments, 
as explained above. 

Figure [5a| shows the success probability, as a func¬ 
tion of temperature, for the D-Wave Two chip, open 
system quantum simulation and the classical-path 
model (SVMC) for = 0.44. The D-Wave Two 
experimental data clearly shows thermal reduction: 
decreasing probability of success (ground state pop¬ 
ulation) with temperature. This is a consequence of 
quantum tunneling, as seen in open quantum sys¬ 
tem theory. Contrary to the experimental data, 
SVMC shows thermal activation, with increasing suc¬ 
cess probability for increasing temperature. For a 


wide range of plausible parameters, only the quantum 
models show thermal reduction in these instances. 
The probability of success of SVMC is also lower than 
D-Wave Two data at the same temperature. 

For /il close to the degeneracy value = J/2 the 
minimum gap becomes small, as seen in Fig. [3] 
Where flio -C W, the adiabatic basis of the instan¬ 
taneous mutiqubit states {|'0o(s)}> lt/h(s))} loses its 
physical significance. Because the coupling to the 
bath is relatively strong here, the system quickly ap¬ 
proaches the states corresponding to predominantly 
opposite cluster orientations, similar to diabatic states 
(see inset of Fig. |3j). Transitions between these states, 
also called pointer states m, occur at a much slower 
rate as a consequence of the polaronic effect. As a 
result, for sufficiently small mininum gaps the mul¬ 
tiqubit freezing starts before the avoided crossing 
and the success probability increases with tempera¬ 
ture [15] . 

A generalization of the 16 qubit problem to a larger 
number of qubits is achieved by studying problems 
that contain the same “motif” multiple times within 
the connectivity graph, see Fig. [6] The success proba¬ 
bilities for up to 200 qubits are shown in Fig. [7] We fit 
the average success probability as p{n q ) oc exp(—cm g ), 
where n q is the number of qubits. The fitting expo¬ 
nent a for the D-Wave Two data is (1.1 ±0.05) • 10 —2 , 
while the fitting exponent for the SVMC numerics is 
(2.8 ± 0.17) • 10~ 2 . We conclude that, for instances 
with multiqubit quantum tunneling, the D-Wave Two 
processor returns the solution that minimizes the en¬ 
ergy with consistently higher probability than physi¬ 
cally plausible models of the hardware that only em¬ 
ploy product states and do not allow for multiqubit 
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(a) Probability vs. temperature. 


(b) Probability vs. /il- 


FIG. 5. Plots (a) and (b) show the probability of success for D-Wave (blue points) versus /il and temperature. The 
decreasing probability with increasing temperature using hi, = 0.44 is only matched with theories based on quantum 
tunneling. This is the opposite tendency to thermal activation (SVMC). Both Redfield and NIB A use only measured 
parameters (no fitting). In this temperature range the lowest two states (the double well potential) account for all the 
probability (0.9998 for D-Wave, 0.99998 for SVMC). 


tunneling transitions. 


The correlation between D-Wave’s experimental 
data and Path Integral Monte Carlo along the Quan¬ 
tum Annealing schedule (PIMC-QA) has been stud¬ 
ied in recent works nn na ns. For complete¬ 
ness, we study PIMC-QA using similar parameters as 
in HS3- PIMC-QA gives a probability of success be¬ 
tween SVMC and the quantum models (Figs. 5b and 
[7]), and does not show thermal reduction for h l = 0.44 
(Fig. [5a|. 


A way to think of multiqubit tunneling as a com¬ 
putational resource is to regard it as a form of large 
neighborhood search. Collective tunneling transitions 
involving K qubits explore a K variable neighbor¬ 
hood, and there is a combinatorial number of such 
neighborhoods. We find that the current generation 
D-Wave Two annealer enables tunneling transitions 
involving at least 8 qubits. It will be an important fu¬ 
ture task to determine the maximal K attainable by 
current technology and how large it can be made in 
next generations. The larger I \, the easier it should 
be to translate the quantum resource “AT-qubit tun¬ 
neling” into a possible computational speedup. We 
want to emphasize that this paper does not claim to 


have established a quantum speedup. To this end one 
would have to demonstrate that no known classical al¬ 
gorithm finds the optimal solution as fast as the quan¬ 
tum process. To establish such an advantage it will 
be important to study to what degree collective tun¬ 
neling can be emulated in classical algorithms such as 
Quantum Monte Carlo or by employing cluster up¬ 
date methods. However, the collective tunneling phe¬ 
nomena demonstrated here present an important step 
towards what we would like to call a physical speedup: 
a speedup relative to a hypothetical version of the 
hardware operated under the laws of classical physics. 
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